%%%%%%%%%%%
%% Forecast ...
%%%%%%%%%%%
close all;
clear all;

%% Read the model
[m,p,mss] = readmodel(false);

%% Load the the historical data from the file '.csv'
% (the file 'kalm_his.csv' is saved by the 'kalmanfilter' program)
% h = dbload('history.csv');
h = dbload('kalm_his.csv');

%% Define the time frame of the forecast
% Change the time frame depending on you data and forecast period!
startfcast = qq(2014,2);
endfcast = qq(2020,4);
fcastrange = startfcast:endfcast;

%% In the forecast simulation variables can be set as numbers (for instance
% for external development as opposed to automatically used AR(1) processes). 
% These preset variables are in the historical database and they require 
% creation of a 'plan' for this variable. The command 'plan' creates plan 
% for model 'm' over the specified horizon in the object 'simplan'. These
% variables must be marked using the command 'exogenize', while the residual
% used for equalizing the model equation for the fixed variable must be marked using
% command 'endogenize'. Note that any variable that is fixed must have residual term
% as a part of its equation. The forecast is simulated in a similar way as 
% without the fixing. See example below where external inflation 
% and external interest rates could be taken, say, from IMF's WEO or 
% Consensus Forecast. Be aware that you must specify exactly till when
% the data is available in the database (line 35).

% xrange = startfcast:qq(2012,4); 
% simplan = plan(m,xrange);
% simplan = exogenize(simplan,{'dot_x_cpi','x_rn'},xrange);
% simplan = endogenize(simplan,{'e_dot_x_cpi','e_x_rn'},xrange);

%% Simulate
% Command 'simulate' simulates the model 'm' based on the database 'h' over
% the forecast range 'fcastrange'. Results are written into the object 's'.
s = simulate(m,h,fcastrange,'nonlinearise=',30,'maxiter=',1000);

%% In case you use exact numbers for external variables, i.e. foreing inflation,
%% interest rate or oil prices, the object "simplan" created below has to be inserted
%% into the command simulate. Comment line 44; uncomment line 50 and the line s = [...] after NTF.

% simplan = plan(m,startfcast:endfcast);

% Specify the exogenize and endogenize commands for each variable, 
% paying attention to each variable periods

%% External variables ("the rest of the world")

% h.dot_x_cpi(startfcast:startfcast+8) = [1,1,1,1,1,1,1,1,1]
% h.x_rn(startfcast:startfcast+3) = [0.25,0.25,0.25,0.25];
% h.lx_gdp_gap(startfcast:startfcast+8) = [-2.7,-2.6,-2.5,-2,-1.5,-1.5,-1,-1,0];
% simplan = exogenize(simplan,{'dot_x_cpi'},startfcast:startfcast+8);
% simplan = exogenize(simplan,{'x_rn'},startfcast:startfcast+3);
% simplan = exogenize(simplan,{'lx_gdp_gap'},startfcast:startfcast+8);
% simplan = endogenize(simplan,{'e_dot_x_cpi'},startfcast:startfcast+8);
% simplan = endogenize(simplan,{'e_x_rn'},startfcast:startfcast+3);
% simplan = endogenize(simplan,{'e_lx_gdp_gap'},startfcast:startfcast+8)

%% Domestic policy variables (e.g., temporary ER peg, "forward guidance")
%% Use sparringly!

% h.ls(startfcast:startfcast+4) = [329,329,329,329,329];
% h.rn(startfcast:startfcast+4) = [0.46,0.46,0.46,0.46,0.46];
% simplan = exogenize(simplan,{'ls'},startfcast:startfcast+4);
% simplan = exogenize(simplan,{'rn'},startfcast:startfcast+4);
% simplan = endogenize(simplan,{'e_ls'},startfcast:startfcast+4)
% simplan = endogenize(simplan,{'e_rn'},startfcast:startfcast+4)

%% Near-term forecast (NTF)

% h.dot_cpi_x(startfcast:startfcast+1) = [0,0];
% simplan = exogenize(simplan,{'dot_cpi_x'},startfcast:startfcast+1);
% simplan = endogenize(simplan,{'e_dot_cpi_x'},startfcast:startfcast+1);
% 
% h.dot_cpi_food(startfcast:startfcast+1) = [0,0];
% simplan = exogenize(simplan,{'dot_cpi_food'},startfcast:startfcast+1);
% simplan = endogenize(simplan,{'e_dot_cpi_food'},startfcast:startfcast+1);
% 
% h.dot_cpi_oil(startfcast:startfcast+1) = [0,0];
% simplan = exogenize(simplan,{'dot_cpi_oil'},startfcast:startfcast+1);
% simplan = endogenize(simplan,{'e_dot_cpi_oil'},startfcast:startfcast+1);

% s = simulate(m,h,fcastrange,'plan',simplan,'nonlinearise=',30,'maxiter=',1000);

%% Uncomment the next line if using a simplan:
% s = simulate(m,h,fcastrange,'plan=',simplan,'nonlinearise=',30,'maxiter=',1000);

% Command 'dbextend' puts together the historical database 'h' and the
% results of the simulation saved in object 's'. Single database 's' is
% created
s.rn = s.rs;
h = dbextend(h,s);

%% Additional forecast-implied data to be saved
% Oil prices, level, USD/barrel
% USD/EUR exchange rate, level
for i = startfcast:endfcast
h.lwoil(i) = h.lwoil(i-1) + h.dot_woil(i)/4;
h.woil(i) = exp(h.lwoil(i)/100);
h.ls_cross(i) = h.ls_cross(i-1) + h.dot_s_cross(i)/4;
h.s_cross(i) = exp(h.ls_cross(i)/100);
end

% Results are saved in file 'fcastdata.mat'  
dbsave('fcastdata.csv',h);

%% Graphs and Tables
% Prepares the forecast report: graphs and tables in the Acrobat .pdf format
Tablerng = startfcast-3:startfcast+7;
Plotrng = startfcast-3:startfcast+11;
Histrng = startfcast-3:startfcast-1;

% Specify country and units for exchange rate
country = 'The Czech Republic';
exchange = 'CZK/EUR';

% Report
x = report.new(country,'visible',true);

% Figures
sty = struct();
sty.line.linewidth = 1.5;
sty.line.linestyle = {'-';'--';':'};
sty.line.color = {'k';'k';'k'};
sty.axes.box = 'on';
sty.legend.location = 'southOutside';
sty.legend.orientation = 'horizontal';

x.figure('Forecats - Main Indicators','subplot',[3,2],'style',sty,'range',Plotrng,'dateformat','YYYY:P');

x.graph('Inflation','legend',true);
x.series('',[h.dot_cpi h.dot4_cpi h.target],'legendEntry=',{'q-o-q','y-o-y','Target'});
x.highlight('',Histrng);

x.graph('Nominal Interest Rate','legend',false);
x.series('',[h.rn]);
x.highlight('',Histrng);

x.graph(['Nominal Exchange Rate - ' exchange],'legend',false);
x.series('',[exp(h.ls/100)]);
x.highlight('',Histrng);

x.graph('Nominal Exchange Rate','legend',true);
x.series('',[h.dot_s (h.ls - h.ls{-4})],'legendEntry=',{'q-o-q','y-o-y'});
x.highlight('',Histrng);

x.graph('Output Gap','legend',false);
x.series('',[h.lgdp_gap]);
x.highlight('',Histrng);

x.graph('Monetary Conditions','legend',true);
x.series('', [h.mci h.rr_gap h.lz_gap],'legendEntry=',{'MCI','RIR gap', 'RER gap'});
x.highlight('',Histrng);

x.pagebreak();

% Tables
TableOptions = {'range',Tablerng,'vline',startfcast-1,'decimal',1,'dateformat','YYYY:P',...
    'long',true,'longfoot','---continued','longfootposition','right'};

x.table('Forecast - Main Indicators',TableOptions{:});

x.subheading('');
  x.series('CPI ',h.dot4_cpi,'units','% (y-o-y)');
  x.series('',h.dot_cpi,'units','% (q-o-q)');
  x.series('Target',h.target,'units','%');
x.subheading('');  
  x.series('Exchange Rate',exp(h.ls/100),'units',exchange);
  x.series('',(h.ls-h.ls{-4}),'units','% (y-o-y)');
x.subheading('');
  x.series('GDP',h.dot4_gdp,'units','% (y-o-y)');
x.subheading('');
  x.series('Interest Rate',h.rn,'units','% p.a.');

x.subheading('');
x.subheading('Inflations');
  x.series('CPI ',h.dot4_cpi,'units','% (y-o-y)');
  x.series('Core Inflation',h.dot4_cpi_x,'units','% (y-o-y)');
  x.series('Food Prices',h.dot4_cpi_food,'units','% (y-o-y)');
  x.series('Oil Prices',h.dot4_cpi_oil,'units','% (y-o-y)');

x.subheading('');
x.subheading('Real Economy');
  x.series('Output Gap',h.lgdp_gap,'units','%');
  x.series('GDP',h.dot_gdp,'units','% (q-o-q)');
  x.series('Potential GDP',h.dot_gdp_eq,'units','% (q-o-q)');
  
x.subheading('');
x.subheading('Monetary Conditions');
  x.series('Monetary Conditions',h.mci,'units','%');
  x.series('Real Interest Rate Gap',h.rr_gap,'units','p.p.');
  x.series('Credit Premium',h.cr_prem,'units','p.p.');
  x.series('Real Exchange Rate Gap',h.lz_gap,'units','%');

x.pagebreak();
x.table('Forecast - Inflation Decomposition',TableOptions{:});

x.subheading('Contributions');
  x.series('CPI',h.dot_cpi,'units','% (q-o-q)');
  x.series('Core Inflation',(1-h.lw_food-h.lw_oil)*h.dot4_cpi_x,'units','p.p.');
  x.series('Food Prices',h.lw_food*h.dot_cpi_food,'units','p.p.');
  x.series('Oil Prices',h.lw_oil*h.dot_cpi_oil,'units','p.p.');
  
x.subheading('');
x.subheading('Core Inflation');
  x.series('Core Inflation',h.dot_cpi_x,'units','%');
  x.series('Lag',h.b1*h.dot_cpi_x{-1},'units','p.p.');
  x.series('Expectations',(1-h.b1)*h.dot_cpi{+1},'units','p.p.');
  x.series('RMC',h.b2*h.rmc,'units','p.p.');
  x.series('RMC - Domestic',h.b2*h.b3*h.lgdp_gap,'units','p.p.');
  x.series('RMC - Imported',h.b2*(1-h.b3)*h.lz_gap,'units','p.p.');
  x.series('Residual',h.e_dot_cpi_x,'units','p.p.');

x.subheading('');
x.subheading('Food Prices');
  x.series('Food Prices',h.dot_cpi_food,'units','%');
  x.series('Lag',h.b21*h.dot_cpi_food{-1},'units','p.p.');
  x.series('Expectations',(1-h.b21)*h.dot_cpi{+1},'units','p.p.');
  x.series('Relative Food Prices',h.b22*h.b23*h.lz_food_gap,'units','p.p.');
  x.series('Business Cycle',h.b22*(1-h.b23)*h.lgdp_gap,'units','p.p.');
  x.series('Residual',h.e_dot_cpi_food,'units','p.p.');
  
x.subheading('');
x.subheading('Oil Prices');
  x.series('Oil Prices',h.dot_cpi_oil,'units','%');
  x.series('Lag',h.b31*h.dot_cpi_oil{-1},'units','p.p.');
  x.series('Expectations',(1-h.b31-h.b32)*h.dot_cpi{+1},'units','p.p.');
  x.series('World Oil Price',h.b32*h.dot_woil,'units','p.p.');
  x.series('Exchange Rate',h.b32*h.dot_s,'units','p.p.');
  x.series('Cross Exchange Rate',h.b32*h.dot_s_cross,'units','p.p.');
  x.series('Residual',h.e_dot_cpi_oil,'units','p.p.');

x.pagebreak();
x.table('Forecast - Demand and Supply',TableOptions{:});

x.subheading('Ouptut Gap Decomposition');
  x.series('Output Gasp',h.lgdp_gap,'units','%');
  x.series('Lag',h.a1*h.lgdp_gap{-1},'units','p.p.');
  x.series('Monetary Conditions',-h.a2*h.mci,'units','p.p.');
  x.series('Real Interest Rate',-h.a2*h.a4*h.rr_gap,'units','p.p.');
  x.series('Credit Premium',-h.a2*h.a4*h.cr_prem,'units','p.p.');
  x.series('Real Exchange Rate',-h.a2*(1-h.a4)*(-h.lz_gap),'units','p.p.');
  x.series('Foreign Output Gap',h.a3*h.lx_gdp_gap,'units','p.p.');
  x.series('Residual',h.e_lgdp_gap,'units','p.p.');
  
x.subheading('');
x.subheading('Supply Side Assumptions');
  x.series('Potential Output',h.dot_gdp_eq,'units','% (q-o-q)');
  x.series('',(h.lgdp_eq-h.lgdp_eq{-4}),'units','% (y-o-y)');
  x.subheading('');
  x.series('Eq. Real Interest Rate',h.rr_eq,'units','%');
  x.subheading('');
  x.series('Eq. Real Exchange Rate',h.dot_z_eq,'units','% (q-o-q)');
  x.series('',(h.lz_eq-h.lz_eq{-4}),'units','% (y-o-y)'); 
  
x.pagebreak();
x.table('Forecast - Policy Decomposition',TableOptions{:});

x.subheading('Interest Rate Decomposition');
  x.series('Interest Rate',h.rn,'units','% p.a.');
  x.series('Lag',h.f1*h.rn{-1},'units','p.p.');
  x.series('Neutral Rate',(1-h.f1)*h.rn_neutral,'units','p.p.');
  x.series('Expected Inflation',(1-h.f1)*h.f2*(h.dot4_cpi{+4} - h.target{+4}),'units','p.p.');
  x.series('Output Gap',(1-h.f1)*h.f3*h.lgdp_gap,'units','p.p.');
  x.series('Residual',h.e_rn,'units','p.p.');

x.subheading('');
x.subheading('Monetary Conditions Decomposition');
  x.series('Monetary Conditions',h.mci,'units','%');
  x.series('Real Interest Rate Gap',h.a4*h.rr_gap,'units','p.p.');
  x.series('Credit Premium',h.a4*h.cr_prem,'units','p.p.');
  x.series('Real Exchange Rate Gap',(1-h.a4)*(-h.lz_gap),'units','p.p');
  
x.pagebreak();
x.table('Forecast - Foreign Variables',TableOptions{:});

x.subheading('European Monetary Union');
  x.series('Inflation',h.dot_x_cpi,'units','% (q-o-q)');
  x.series('Interest Rate',h.x_rn,'units','% p.a.');
  x.series('Output Gap',h.lx_gdp_gap,'units','%');

x.subheading('');
x.subheading('World Food and Oil Prices');
  x.series('IMF Food Price Index',h.dot_wfood,'units','% (q-o-q)');
  x.series('Oil Prices',h.dot_woil,'units','% (q-o-q)');
  x.series('',h.woil,'units','USD/Barrel');

x.subheading('');
x.subheading('USD/EUR Exchange Rate');
  x.series('',h.dot_s_cross,'units','% (q-o-q)');
  x.series('',h.s_cross,'units','level');

x.publish('Forecast','display',false);
disp('Done!');
    
